In [1]:

    

import matplotlib.patches as mpatches
import matplotlib.pyplot as plt
import shapely.geometry as sgeom
import cartopy
import cartopy.crs as ccrs
import cartopy.io.shapereader as shpreader


    

In [2]:

    

figure(figsize=(8,6))
ax = plt.axes(projection=cartopy.crs.PlateCarree())

ax.add_feature(cartopy.feature.LAND)
ax.add_feature(cartopy.feature.OCEAN)
ax.add_feature(cartopy.feature.COASTLINE)
ax.add_feature(cartopy.feature.BORDERS, linestyle=':')
ax.add_feature(cartopy.feature.LAKES, alpha=0.5)
ax.add_feature(cartopy.feature.RIVERS)

ax.set_extent([-20, 60, -40, 40])
title('Africa in PlateCarree Projection')
plt.show()


    

In [3]:

    

def waves_sample_data(shape=(73, 145)):
    """Returns ``lons``, ``lats`` and ``data`` of some fake data."""
    nlats, nlons = shape
    lats = np.linspace(-np.pi / 2, np.pi / 2, nlats)
    lons = np.linspace(0, 2 * np.pi, nlons)
    lons, lats = np.meshgrid(lons, lats)
    wave = 0.75 * (np.sin(2 * lats) ** 8) * np.cos(4 * lons)
    mean = 0.5 * np.cos(2 * lats) * ((np.sin(2 * lats)) ** 2 + 2)

    lats = np.rad2deg(lats)
    lons = np.rad2deg(lons)
    data = wave + mean

    return lons, lats, data


    

In [4]:

    

figure(figsize=(8,4))
ax = plt.axes(projection=ccrs.Mollweide())

lons, lats, data = waves_sample_data()

ax.contourf(lons, lats, data,
            transform=ccrs.PlateCarree(),
            cmap='spectral')
ax.coastlines()
ax.set_global()
title('Planetary Waves on Mollweide Projection')
plt.show()


    

In [5]:

    

def katrina_sample_data():
    """
    Returns a list of latitudes and a list of longitudes (lons, lats)
    for Hurricane Katrina (2005).

    The data was originally sourced from the HURDAT2 dataset from AOML/NOAA:
    http://www.aoml.noaa.gov/hrd/hurdat/newhurdat-all.html on 14th Dec 2012.

    """
    lons = [-75.1, -75.7, -76.2, -76.5, -76.9, -77.7, -78.4, -79.0,
            -79.6, -80.1, -80.3, -81.3, -82.0, -82.6, -83.3, -84.0,
            -84.7, -85.3, -85.9, -86.7, -87.7, -88.6, -89.2, -89.6,
            -89.6, -89.6, -89.6, -89.6, -89.1, -88.6, -88.0, -87.0,
            -85.3, -82.9]

    lats = [23.1, 23.4, 23.8, 24.5, 25.4, 26.0, 26.1, 26.2, 26.2, 26.0,
            25.9, 25.4, 25.1, 24.9, 24.6, 24.4, 24.4, 24.5, 24.8, 25.2,
            25.7, 26.3, 27.2, 28.2, 29.3, 29.5, 30.2, 31.1, 32.6, 34.1,
            35.6, 37.0, 38.6, 40.1]

    return lons, lats


    

In [6]:

    

figure(figsize=(8,6))
ax = plt.axes([0.01, 0.01, 0.98, 0.98], projection=ccrs.PlateCarree())

ax.set_xlim([-125, -66.5])
ax.set_ylim([20, 50])

states_shp = shpreader.natural_earth(resolution='110m',
                                     category='cultural',
                                     name='admin_1_states_provinces_shp')

lons, lats = katrina_sample_data()

# to get the effect of having just the states without a map "background"
# turn off the outline and background patches
ax.background_patch.set_visible(False)
ax.outline_patch.set_visible(False)

plt.title('US States which intersect the track '
          'of Hurricane Katrina (2005)')

# turn the lons and lats into a shapely LineString
track = sgeom.LineString(zip(lons, lats))

# buffer the linestring by two degrees (note: this is a non-physical
# distance)
track_buffer = track.buffer(2)

for state in shpreader.Reader(states_shp).geometries():
    # pick a default color for the land with a black outline,
    # this will change if the storm intersects with our track
    facecolor = [0.9375, 0.9375, 0.859375]
    edgecolor = 'black'

    if state.intersects(track):
        facecolor = 'red'
    elif state.intersects(track_buffer):
        facecolor = '#FF7E00'

    ax.add_geometries([state], ccrs.PlateCarree(),
                      facecolor=facecolor, edgecolor=edgecolor)

ax.add_geometries([track_buffer], ccrs.PlateCarree(),
                  facecolor='#C8A2C8', alpha=0.5)
ax.add_geometries([track], ccrs.PlateCarree(),
                  facecolor='none')

# make two proxy artists to add to a legend
direct_hit = mpatches.Rectangle((0, 0), 1, 1, facecolor="red")
within_2_deg = mpatches.Rectangle((0, 0), 1, 1, facecolor="#FF7E00")
labels = ['State directly intersects\nwith track',
          'State is within \n2 degrees of track']
plt.legend([direct_hit, within_2_deg], labels,
           loc=3, fancybox=True)

plt.show()


    

In [7]:

    

figure(figsize=(12,12))
rotated_pole = ccrs.RotatedPole(pole_latitude=45, pole_longitude=180)

box_top = 45
x, y = [-44, -44, 45, 45, -44], [-45, box_top, box_top, -45, -45]

ax = plt.subplot(211, projection=rotated_pole)
ax.stock_img()
ax.coastlines()
ax.plot(x, y, marker='o', transform=rotated_pole)
ax.fill(x, y, color='coral', transform=rotated_pole, alpha=0.4)
ax.gridlines()
title('Rotated Pole Coordinates -- Tough in Basemap')


ax = plt.subplot(212, projection=ccrs.PlateCarree())
ax.stock_img()
ax.coastlines()
ax.plot(x, y, marker='o', transform=rotated_pole)
ax.fill(x, y, transform=rotated_pole, color='coral', alpha=0.4)
ax.gridlines()
plt.show()


    

In [7]: